<!DOCTYPE group PUBLIC "-//W3C//DTD XHTML 1.0 Transitional//EN"
"http://www.w3.org/TR/xhtml1/DTD/xhtml1-transitional.dtd">
<group>

%tableofcontents;

<p>This tutorial shows how to estiamte <a href="%dox:gmm;">Gaussian
mixture model</a> using the VlFeat implementation of
the <em>Expectation Maximization</em> (EM) algorithm.</p>

<p>A GMM is a collection of $K$ Gaussian distribution. Each
distribution is called a <em>mode</em> of the GMM and represents a
cluster of data points.  In computer vision applications, GMM are
often used to model <em>dictionaries of visual words</em>. One
important application is the computation
of <a href="%pathto:tut.fisher;">Fisher vectors encodings</a>.</p>

<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->
<h1 id="tut.gmm.introduction">Learning a GMM with expectation maximization</h1>
<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->

<p>Consider a dataset containing 1000 randomly sampled 2D points:</p>

<precode type='matlab'>
numPoints = 1000 ;
dimension = 2 ;
data = rand(dimension,N) ;
</precode>

<p>The goal is to fit a GMM to this data. This can be obtained by
running the <code>vl_gmm</code> function, implementing
the <a href="%dox:gmm-em;">EM algorithm</a>.</p>

<precode type='matlab'>
numClusters = 30 ;
[means, covariances, priors] = vl_gmm(data, numClusters) ;
</precode>

<p>Here <code>means</code>, <code>covariances</code>
and <code>priors</code> are respectively the means $\mu_k$, diagonal
covariance matrices $\Sigma_k$, and prior probabilities $\pi_k$ of
the <code>numClusters</code> Gaussian modes.</p>

<p>These modes can be visualized on the 2D plane by plotting ellipses
corresponding to the equation:
\[
   \{ \bx: (\bx-\mu_k)^\top \Sigma_k^{-1} (\bx-\mu_k) = 1 \}
\]
for each of the modes. To this end, we can use
the <code>vl_plotframe</code>:</p>

<precode type='matlab'>
figure ;
hold on ;
plot(data(1,:),data(2,:),'r.') ;
for i=1:numClusters
    vl_plotframe([means(:,i)' sigmas(1,i) 0 sigmas(2,i)]);
end
</precode>

<p>This results in the figure:</p>

<div class="figure">
  <image src="%pathto:root;demo/gmm_2d_rand.jpg"/>
  <div class="caption">GMM fittting 2D random points.</div>
</div>

<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->
<h2 id="tut.gmm.cov">Diagonal covariance restriction</h2>
<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->

<p>Note that the ellipses in the previous example are axis
alligned. This is a restriction of the <code>vl_gmm</code>
implementation that imposes covariance matrices to be diagonal.</p>

<p>This is suitable for most computer vision applications, where
estimating a full covariance matrix would be prohebitive due to the
relative high dimensionality of the data. For example, when clustering
SIFT features, the data has dimension 128, and each full covariance
matrix would contain more than 8k parameters.</p>

<p>For this reason, it is sometimes desirable to globally decorrelated
the data before learning a GMM mode. This can be obtained by
pre-multiplying the data by the inverse of a square root of its
covariance.</p>

<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->
<h1 id="tut.gmm.initialization">Initializing a GMM model before running EM</h1>
<!-- ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ -->

<p>The EM algorithm is a local optimization method, and hence
particularly sensitive to the initialization of the model. The
simplest way to initiate the GMM is to pick <code>numClusters</code>
data points at random as mode means, initialize the individual
covariances as the covariance of the data, and assign equa prior
probabilities to the modes. This is the default initialization
method used by <code>vl_gmm</code>.</p>

<p>Alternatively, a user can specifiy manually the initial paramters
of the GMM model by using the <code>custom</code> initalization
method. To do so, set
the <code>'Initialization'</code> option  to <code>'Custom'</code> and
also the options <code>'InitMeans'</code>, <code>'InitCovariances'</code> and
<code>'IniPriors'</code> to the desired values.</p>

<p>A common approach to obtain an initial value for these parameters
is to run KMeans first, as demonstrated in the following code
snippet:</p>

<precode type='matlab'>
numClusters = 30;
numData = 1000;
dimension = 2;
data = rand(dimension,numData);

% Run KMeans to pre-cluster the data
[initMeans, assignments] = vl_kmeans(data, numClusters, ...
    'Algorithm','Lloyd', ...
    'MaxNumIterations',5);

initCovariances = zeros(dimension,numClusters);
initPriors = zeros(1,numClusters);

% Find the initial means, covariances and priors
for i=1:numClusters
    data_k = data(:,assignments==i);
    initPriors(i) = size(data_k,2) / numClusters;

    if size(data_k,1) == 0 || size(data_k,2) == 0
        initCovariances(:,i) = diag(cov(data'));
    else
        initCovariances(:,i) = diag(cov(data_k'));
    end
end

% Run EM starting from the given parameters
[means,covariances,priors,ll,posteriors] = vl_gmm(data, numClusters, ...
    'initialization','custom', ...
    'InitMeans',initMeans, ...
    'InitCovariances',initCovariances, ...
    'InitPriors',initPriors);
</precode>

<p>The demo scripts <code>vl_demo_gmm_2d</code>
and <code>vl_demo_gmm_3d</code> also produce cute colorized figures
such as these: </p>

<div class="figure">
  <image src="%pathto:root;demo/gmm_2d_shell.jpg"/>
  <div class="caption">The figure shows how the estimated gaussian
  mixture looks like with and without the kmeans initialization.</div>
</div>

</group>

